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ABSTRACT 

Context. Radial velocity surveys such as the Radial Velocity Experiment (RAVE) provide us with measurements of hundreds of 
thousands of nearby stars most of which belong to the Galactic thin, thick disk or halo. Ideally, to study the Galactic disks (both thin 
and thick) one should make use of the multi-dimensional phase-space and the whole pattern of chemical abundances of their stellar 
populations. 

Aims. In this paper, with the aid of the RAVE Survey, we study the thin and thick disks of the Milky Way, focusing on the latter. We 
present a technique to disentangle the stellar content of the two disks based on the kinematics and other stellar parameters such as the 
surface gravity of the stars. Using the Padova Galaxy Model, we checked the ability of our method to correctly isolate the thick disk 
component from the Galaxy mixture of stellar populations. 

Methods. We introduce selection criteria in order to clean the observed radial velocities from the Galactic differential rotation and to 
take into account the partial sky coverage of RAVE. We developed a numerical technique to statistically disentangle thin and thick 
disks from their mixture. 

Results. We deduce the components of the solar motion relative to the Local Standard of Rest (LSR) in the radial and vertical direction, 
(vjy, v w ) Q = (9.87 ± 0.37, 8.01 ± 0.29) km s _1 , the rotational lag of the thick disk component relative to the LSR v thick i ag = 49 km s~ , 
and the square root of the absolute value of the velocity dispersion tensor for the thick disk alone: cr RR = (56.1 ± 3.8)km s , 
a Rir = (29.4±17.2)km s"',cr fc = (10.1±3.3)km s" 1 , <r w = (46.1±6.7)km s -1 ,^ = (5.8±5.1)km s -1 , cr H = (35.1±3.4)km s _1 .The 
analysis of the thin disk is presented in another paper. We find good agreement with previous independent parameter determinations. 
In our analysis we used photometrically determined distances. In the Appendix we show that similar values can be found for the thick 
disk alone as derived in the main sections of our paper even without the knowledge of photometric distances. 

Key words. Stellar kinematics and dynamics - Methods: analytical, numerical - Surveys - Stars: kinematics - Galaxy: structure and 
evolution, thick disk 



1. Introduction 

Stellar radial veloc ity surveys such as t he Radial Velocity 
Experiment RAVE dSteinmetz et all l2006t IZ witter et all 120081: 
Sieb ert et al.l201 lbh can provide, for aU the observed stars, three 
of the six components of the Galactic stellar phase space - two 
directional components and one velocity component. The lack 
of the other three components makes it difficult to fully char- 
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acterize the dynamical behaviour of our Galaxy and its present 
phase-space description from such surveys. 

Despite these limitations, RAVE has the unique merit of 
systematically measuring stellar radial velocities with a pre- 
cision never previously realized for such a large sample of 
stars, providing us with a statistically significant amount of data 
that continues to grow, recently surpassing the 500,000 spec- 
tra. RAVE is a southern-hemisphere survey obtaining spectra 
in the near-infrared Call triplet region of magnitude selected 
stars (9 < / < 13). The usefulness of this vast dataset has 
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already been demonstrated in the context of kinematics in or- 
der to deduce t he characteristics of the stellar velocity dis- 
tributions (e.g.. | VeItz et all 120081: ICas etti-Dinescu et al. 201 1; 
Siebe rt et al l l2008t lHahn et al] 1201 it ICoskunoglu et alJboi ll 
Sieb ert et al] 1201 lap | Wilson et al] 1201 U iKaratas & Klementl 
20121: iBiliret al] 120 12l)T to identify the stellar streams (e.g. , 

aLl boOSllWilhams et alll2011t 



Seabroke et al. 2008; Klement et ; 

Kiss et all 1201 1[). and to study the chemistry of MW com- 



ponents (elg^ iRuchti et al]|2010l 1201 lbt IFu lbright et all 120101: 



lBoecheetal.ll201lt Coskunoglu et al .1 120121)" Here, we carry out 
an ana lysis on topics recently addresse d by ICoskunoglu et al] 
d2011l) and ICasetti-Dinescu et all d2011l) . leaving aside for the 
moment the int erpretative tools ba sed on the orbit integration as 
developed by Wilson et al] (1201 ll) . In other words, although we 
will make use of a Galactic potential model, the evolution with 
time t will not be the direct object of our a nalysis. The mag nitude 
range of RAVE (9 < / < 13) (see e.g ., ISteinmetz etal]|2006t 
Zwit ter et al.ll2008t Isiebert et al]|201 lal) implies that most of the 
stars targeted by RAVE belong to the thin or the thick disk of 
the MW. Moreover, working with RAVE data to investigate the 
kinematics of our Galaxy in a statistical sense, basically means 
deriving the time-independent single-component disk-like dis- 
tribution function (DF) / (x, v), of the stars in the phase space 
mapped with coordinates {x, v), or an overlapping set of DFs 
/ tot (x, v) = 2 /; (x, v) for each stellar population /' that we are 
able to disentangle. In terms of direct observations of the phase 
space, and within reasonable errors, RAVE provides the direc- 
tion (line of sight, l.o.s.) and the velocity along the l.o.s. for each 
star. The radial velocity is the observed projection of the true 
heliocentric observed velocity vector Vhd of the star along the 
line of sight, v r = (Vh e i,?h e i) (with (..., ...} denoting the standard 
inner product). The direction f h e i = X - x /||x - x Q || is just the 
heliocentric position vector divided by the unknown heliocentric 
distance of each star, ||rh e il| = ||x - x Q || = rh e i (with components 
{x,y,z} in a suitable reference system, see Section©. 

The analysis of the kinematic properties of the thin and thick 
disks is of paramount importance for the comprehension of the 
origin and formation of our Galaxy as well as of every disk 
galaxy. The object of our analysis in the present paper is mainly 
the thick disk. A separate paper is dealing w ith the thin disk 
analysis (hereafter Paper II, Pasetto et al1 l2012l) . 

In this paper we investigate the kinematic properties of the 
Galactic disk components by applying a methodology based 
on Singular Value Decomposition (SVD). SVD permits one 
to find the solution, in a least-squares sense, of an inhomo- 
geneous system of N linear equations in the case of singular 
matrices (where N is the number of stars we select). A num- 
ber of studies have used this method in the context of proper 
motion survey analyses: see, e.g., a fully analytical exercise 
de veloped for HIPPARCOS proper motion and p arallax data 



Dehnen & Binnevl dl998l) or a recent work by iFuchs et aD 
(2009) combinin g data release 7 of the Sloan Digital Sky Survey 
dAbazaiian et al] [20091) with astrometric USNO-B data supple- 
mented by photometric distances. 

Finally, several authors devise d methods to derive photo- 
metric distances for RAVE (e.g.. i Klement et all I2008L 1201 1; 
IZwitter et al.1 l2010t iBreddels et ai1l2010t iBurnett et alJl201lh . 
Here we wil l adop t the internal release of the catalogue of 
IZ witter et al.l d2010b which includes different distance deter- 
minations dependin g on different sets of isochrones for about 
260,000 stars (see IZ witter et al.l I2010L for an extended discus- 
si on). We will use the Yon sei-Yale distances based on isochrones 
of lDemarque et al] d2004 . 



The structure of this paper is as follows: in Section|2]we will 
prepare the data for the analysis by cleaning the radial veloci- 
ties from the effects of Galactic rotation. In Section[3]we briefly 
review the inversion techniques we are going to apply for the 
analysis of the data. In Section|4]we present the results and dis- 
cuss them in Section© 

2. Preparing the data 

Our procedure requires the following steps. We need to 

1. correct the RAVE survey data for large scale effects (this 
Section©, 

2. isolate the thick disk (see Section©, 

3. evaluate the thick disk kinematics (see Section©. 

As an additional result of this procedure, we will be able to esti- 
mate two of the three components of the solar motion relative to 
the LSR and the velocity lag of the thick disk component. 

In this way, we can determine some of the moments of 
the underlying multi-component distribution function / tot (x, v). 
While strictly speaking only the infinite series of the moments is 
equivalent to the original distribution function, only the first few 
moments will be the subject of our exercise. Thus, our study is 
intended to be of an exploratory nature rather than exhaustive, 
for a set of problems that are far from being fully solved. 

2.1. Radial velocity component from the Galaxy differential 
rotation 

For a given position of a star in the Galactic disk, we want to 
estimate the influence of the differential rotation of the Galaxy 
in the radial velocity component. This is necessary to minimize 
the influence of the radial component due to the Galactic rotation 
on the observed radial velocity. 

To compute this first step we note that, for a star at a given 
position x in the Milky Way's disk, the following general vector 
relation holds: 



VLSR (Xq) + Vq + Vhel = V c (X) + V„ 



(1) 



where \ c (x) is the mean rotational velocity at the given position 
x in the Galactic disk reference system centred on the bary centre 
of the Milky Way (MW) (considered at rest or in rectilinear un- 
perturbed motion). Vlsr (x©) is the Solar Local Standard of Rest 
which differs from the Local Standard of Rest (LSR) speed at 
any other location in the Galaxy, v© is the Solar peculiar veloc- 
ity relative to the Solar LSR and \ p is the peculiar velocity of the 
star relative to its own mean Galactic rotational speed. 

It is also convenient to introduce here the orthonormal stan- 
dard system of reference in the velocity space (O; U, V, W). This 
reference system is centred on the velocity of the Solar LSR, O, 
with U aligned with the reference system in the configuration 
space pointing to the Galactic centre, V aligned with the rota- 
tion of the Galaxy and W pointing to the north Galactic pole 
(NGP). We will call the generic velocity components in this ref- 
erence system v = {vu, vy, vw}- From Eqn. (©, if we express the 
component of the radial velocity due to Galactic rotation as a 
function of the star's longitude and latitude (I, b), we get 

v? = (llv c («)|| y - Vlsr (Re)) cos b sin /. (2) 

For simplicity we assume cylindrical symmetry {0;R,<p,z) in 
the configuration space for the Galactic model, Vlsr (xq) = 
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{0, Vlsr (Rq) , 0). v = {vu,Q,vv,e,vw,e} are the components of 
the peculiar motion of the Sun relative to the Solar LSR. R is the 
Galactocentric radial distance in cylindrical coordinates for the 
position of a given star. R G is the position of the Sun on the plane 
of symmetry of the MW at an azimuthal position (p Q = 0. The 
height of the Sun relative to the plane of symmetry of the Galaxy 
is hereafter neglected, z = 0. The orthonormal rotation matrix 
R r R = 1, i.e. R e 5(9(3) with det = +1 is adopted here to 
transform a velocity vector (given in terms of radial velocity v r 
and motion along the Galactic coordinates (v/, v^,)) into a vector 
in the velocity space reference system vector (O; U, V, W). This 
orthonormal rotation matrix is defined as 



R = 



cos b cos/ cos/? sin/ sin/?' 

- sin / cos / 

- sin b cos / - sin b sin / cos b, 



(3) 
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2.2. The influence of distance 

We show now how the influence of the Galactic rotation on the 
radial velocity is strongly dependent on the viewing direction, 
i.e., the direction (/, b) for a specific star of the survey, and on 
the mean streaming velocity ||v c (x)|| that will be our free tuning 
parameter. On the other hand, we can choose suitable directions 
(/, b) for which the influence of the Galactic rotation on the radial 
velocity is weakly dependent on the distance rh e i of the star along 
(/, b). This result holds for the range of distances from the Sun in 
which RAVE dwarf stars are mainly located, say r^i < l.Okpc 
as shown in Appendix C. 

In equation (0 the Galactic radial distance of a template star 
can be expressed as a function of the heliocentric distance, H, e i. 
In particular the the first term in the sum of Eqn. (0 reads 

R 

l|v c (x)|| -J^ cos Z? sin/ = 
R 

||v e (x)||cos/?sin/ (4) 
^1 + ^-cosZ?(^cosZ?-2cosQ 

We will treat the terms ||v c (x)|| independently later. For simplic- 
ity, we set e = ^ . For our approximation we only consider the 
kinematics in the solar neighbourhood, defined, e.g., as the zone 
where e = ^ = I. Then it follows from equations (0 and (0 
that in this neighbourhood, by expanding the previous equation 
in a Maclaurin series in e, we obtain 

cos b sin / 
VI + £ cos b(e cos b - 2 cos/) 
cos b sin /+ 

ecos 2 /? cos / sin / + o (s 2 ) 

Thus this expansion suggests that is primarily dominated by 
the direction of observation (/, b) and the v c value (which differs 
from star to star and from the thin disk to the thick disk), but 
we expect no strong dependence on the heliocentric distance of 
the star within our neighbourhood. This is supporting the idea 
that the error affecting the photometrically determined distances 
will not be dramatically relevant for the thick disk parameter 
determinations (see also ICasertano et a"flll990l) . At present, we 
must consider the term v c as a free parameter. We postpone 
its treatment and its dependence on the heliocentric distance 
v c (x) = v e (ivi) until Section |U because its behaviour will be 
treated in the context of the inversion technique. Different argu- 
ments hold for the thin disk analysis (see Paper II) where the 



Fig. 1. Isocontour plot on an Aitoff projection of the differences 
in radial velocity due to Galactic rotation 6vf from Eqn. (0 
for test values of v c calculated using Eqn. (l24l and the values 
in Table [T] The gradients of the v elocity dispersion tensor are 
treated as in Vallenari et al. (2006). The assumed values for the 
motion of the Sun relative to the LSR are from this paper. The 
yellow points show the RAVE data selection from Section 12.11 
The black contour lines represent the limit of 14 km s adopted 
in our paper. 

photometric distance errors are the main source of uncertainties 
for the trend of the velocity ellipsoid in the meridional plane. 

We can test the validity of the previous approximation with- 
out using the expansion into a series, but instead by numerically 
solving for the dependence = (rh e i) at different values 
of (/,/?). In particular, we are interested in controlling the er- 
rors introduced by large scale effects such as differential rotation 
and photometric distances for selected distances. The dispersion 
of the thick disk velocity ellipsoid is generally assumed to be 
constant and thus should not introduce any adverse large scale 
effects on the small volume considered. Nevertheless, we will 
check this assumption "a posteriori" in Section l4.4.2l If we limit 
ourselves to consider the RAVE stars within a sphere of radius 
1.0 kpc centred at the Sun's location, from Eqn.© we can cal- 
culate the value of the radial velocity due to the Galactic rotation 

along the line of sight e [0, 1] kpc for each direction in 
the sky specified by (/, b). To evaluate this equation, we test our 
initial guess valu es for Eqn. (0 with th e help of a Galactic pop- 
ulation model by Vallenari et al. (2006) where the MW Galactic 
potential and velocity dispersion profiles are consistently deter- 
mined (see also our Eqn. (l24l and Section !?. 4. U . 

This value of the Galactic component of the radial veloc- 
ity ranges from a maximum, max(v^), and a minimum value, 

min (vj 3 ), whose difference we will call 

6v^ = Imax v!? - min v?| „ ,, (5) 

r \ i r l, V] e[0,l]. 

This difference 6vf can be evaluated at different values of (/, b) 
and for an optimized value of v c . In this way we obtain the back- 
ground in Figure Q~]that shows the Aitoff sky projection contour 
plot obtained from Eqn. (0. As we can see from this contour 
plot, the highest 6v^ (yellowish zones) occur at low latitude and 
the smallest in the green zones. For b — 0, neglecting the motion 
of the Sun relative to t he Solar LSR, Eqn. d5} is ju st the terminal 
velocity relation (e.g., iBinnev & Merrifieldlll998l) . The highest 
induced error is \6vf\ = 16km s _1 at a distance of l.Okpc. This 
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may seem quite large if compared with the mean radial velocity 
error of RAVE data but this error occurs at low latitude regions 
that are not of interest for RAVE because these are not covered 
by the survey (see Figure [TJ. Also the dependence of \ c on the 
range of interest of ivi is again more relevant at low latitude. 

2.3. RAVE data selection criteria 

It becomes evident from the previous considerations that when 
we introduce an indetermination in the distance d, of the /'''-star 
due to the error in the photometrically determined distance Ac/,, 
we are also introducing an error in the correction for the radial 
velocity contribution due to the Galactic rotation. If we proceed 
with a Monte Carlo approach by randomly generating a dis- 
tance d, within Ad, for the /'''-star, the correction vf, for this star 

will assume only values between the two extremes max (v?) 

r hel e[0,l]V r I 

and min (v?) in Eqn. 01). If this correction is smaller than 

r hc ,e[0,l]V > 

the intrinsic error on the radial velocity as derived by RAVE, 
Av* AVE (/;, bi), i.e., Av^ VE (/,-, b-) > <5vf (/;, b t ), we can either ap- 
ply this correction or we can safely neglect it. This is because 
its contribution to v G is too small to bias our consideration, i.e. 

r,i 

its effect cannot be accounted for within the error AvJ^^ (/,-, bi) 
of the observed stellar radial velocity. Equivalently, by shifting 
the /'''-star by its distance uncertainty, it does not affect the radial 
velocity contribution due to the Galactic rotation because it 
is negligible compared to the intrinsic error Av^ VE (/,■, bi). We 
point out that the observational errors, Av r , for which we use 
capital delta (A), are of a different nature than the errors dvf de- 
scribed by Eqn. (0, which have different systematic trends and 
an origin that can b e accounted for with a complete Galactic 
dynamics model (see IVallenar i et al. 2006, their Appendix), and 
our Eqn. (124T). 

In order to adopt an even more restrictive condition, we will 
use all those stars from the RAVE catalogue that have an intrinsic 
error AvJ* AVE (/,, bi) greater than twice the error introduced by 
neglecting the stars' distance dependence on the velocity due to 
the Galactic rotation. In other words, for every star i, we require: 



Avf AVE (h,bi)> 2.06v^ (li,bi). 



(6) 



A test of the validity of all these arguments to recover the correct 
results has been performed on a mock catalogue as explained in 
Sectiongj] 

From the sample of stars that survive the selection defined 
by Eqns. (0 and © shown in Figure [1] we proceed further with 
a few extra selection criteria. We consider stars with signal-to- 
noise ratio S /N > 20. We require the near-infrared colours of the 
RAVE stars taken from 2MASS observations to be in the range 
of /— K G [0.2, 1 . 1] in order to clean the data set from extremely 
young stars. Furthermore, we require vf AVE < 300kms~' to 
avoid stars dynamically not representative of the thick disk 
within 2 kpc from the Sun's position, Av^ AVE < 14kms _1 
to reduce the propagation error while retaining a considerable 
amount of stars (see the velocity contour highlighted in black in 
Fig-UK \b\ > 10° because we are not interested in thin disk stars, 
and log 10 g > 3.5 to avoid the influence of giant stars which may 
sample a much more distant part of the Galaxy. Together with 
Eqn. (O the last criterion represents the most severe cut on the 
total number of stars that we can use. Where possible a further 
cut [Fe /H] > - 1 dex was applied adopting the metallicity deter- 
mination described by Zwitte ret al.l (120081) . 

Starting with an i nitial number of rou ghly 260,000 stars 
in the data release of IZ witter et al.1 ([2010), the final, remain- 



ing number of stars for which we can perform the analysis is 
N tot - 38, 805 which is orders of magnitude larger than the sam- 
ples currently available in the literature (see Fig.[3]for a plot of 
the error distribution). 



3. Inversion techniques for radial velocities 

In the following Section [3~T1 we will present the technique to de- 
termine the principal moments of a DF for a mixture of thin and 
thick disk stars from a radial velocity survey. With a few extra 
assumptions we will be able to disentangle the two components 
and focus our attention on the thick disk as will be shown in 
SectionEHl 



3. 1 . The mixture of the thin and thick disk distribution of 
moments 

Once we have selected the sample of RAVE stars that satisfy 
the conditions laid out in Section I2.3I (see yellow dots in Fig. 
Q}, we can deduce the first moments of the underlying mixture 
distribution function with simple algebra based on a very popular 
technique for dealing with sets of linear equations that can be 
written in the matrix form Ax = b, for the unknown vector x. In 
the following we recall the basics applied to our specific case of 
a radial velocity survey. For more details we refer to Appendix 
A. 

Based on Eqn. (Q} we can express the components of the 
radial velocity vector, say vy, by 



v r cos b COS A 
v r cos b sin I 
v r sin b 



(7) 



where vy = (v, f ) f = pv with p being the idempotent matrix of 
the projection operator along the line of sight. Proceeding com- 
ponent by component we obtain: 



V|| = pVhd = I"" Vhel 



' x 2 


xy 


xz N 






xy 


f 


yz 




vv 


, xz 


n 









(8) 



where (...)■(...) is the ordinary matrix product, a®" is the standard 
tensor power of the generic vector a (and n e N any non-negative 
integer) and f is again the unitary vector of a star whose com- 
ponents are {x, y, z] in the configuration space collinear with the 
velocity space (O; U, V, W) and centred on O. Clearly the matrix 
p is singular (its detp = 1 - ||r|| =1-1=0) thus not permitting 
us to determine its inverse. 

Nevertheless, we can proceed by taking into consideration 
the whole subsample of the selected RAVE data in order to statis- 
tically invert the overdetermined system Ax = b that we obtain 
by defining A as the block diagonal matrix of all the projection 
operators p once the off-diagonal blocks are small (see Appendix 
A for further details). We call this (3Nx3) block matrix of all the 
projection operators simply P. b will be defined from the vectors 
of the observed radial velocities for the selected sample, which 
we call for simplicity again vy. We obtain in this way the system 
Pvhei = V[| that we multiply as usual by the transposed matrix 
P T in order to obtain the square matrix product P r P that can be 

inverted in order to give Vhe l = (P 7 ^) P^ ll = P +V ll> where P + 
is the pseudo inverse of P (iPenrosd Il956l 119551) . Probably the 
most widely used orthogonal decomposition suitable for the so- 
lution of our system is the Si ngular Value Dec omposition due 
to its numerical stability (e.g. JPress et al.l [l986). which we will 
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also adopt here (e.g., by writing P + = Jv^q with W and U 
as orthogonal ma trices and Q as a diagonal matrix (see e.g., 
iPress et al] dl98 6) for an extended discussion and implementa- 
tion techniques). Finally, in order to check the consistency of 
this approximation we will keep track of the sky distribution by 
weighting the obtained matrix P with the isotropic case laid out 
in Appendix A (see Sect. l4.3l and Eqn. A. 3). 

There is an implicit hypothesis that is underlying our ap- 
proach. When we are considering the system of equations Pvhei, 
we need to implicitly assume that either the group of stars we 
are sampling is so extremely local that it is not (or only weakly) 
influenced by the spatial trend of the velocity dispersion tensor, 
or that the population we are sampling has an isothermal disper- 
sion tensor within the distances sampled. In the latter case the 
sample does not need to be local. For the thick disk stars we can 
exploit th e latter, while the first ha s been already exploited in the 
work of Dehnen & Binnev (Il998h . A bi-dimensional schematic 
representation of this assumption can be seen in Fig. [2] the light 
blue arrows represent the unknown velocity vectors of which we 
observe the radial velocities (green arrows) along a given direc- 
tion (grey arrows). The velocity dispersion described by the light 
blue arrows is found (see below) not to depend on the position 
inside or outside the solar radius R G throughout the volume cov- 
ered by our sample. For comparison see Figure 2 of Paper II 
where the trend of the thin disk velocity dispersion tensor along 
the meridional plane is considered. 

This complex procedure provides our first moment of the 
composite distribution function, hereafter simply \. With this 
approach the solar motion relative to a selected sample of stars 
is simply the mean motion of the stars relative to the Sun with 
a change in sign for the U and W directions, say (vu, Vyp) = 
~(yUiVw)- I n me V direction the lag of the component motion 
has to be taken into account with an extra v\ ag whose contribu- 
tion overlaps with v^e to form the mean value Vy. 

In the same way we can get the second, third, and fourth 
central moment, by defining v' = V — V and averaging, e.g., for 
the moment of order 2: 

= (Pv'f 2 = P®V® 2 . (9) 
Hereafter we exploit the symmetries of the matrices: 

oj = i(P®P) + (P®P) r ■ <T 2 EE M [2] • <T 2 (10) 

where er 2 is the dispersion tensor of the radial velocity and er 2 
is the ordinary pressure tensor. In the same way, we define the 
projected third moment along the l.o.s. as 

S,| ee = pga v /83 = m [3] . S , (11) 

and the projected fourth moment as 

T|, = = P®V»4 = M [4] . T, (12) 

where the suitable symmetries have been exploited. The inter- 
ested reader can find the element by element characterization of 
the matrices M[,j in Appendix A. 

The errors in these quantities are given by standard statisti- 
cal tools coming from linear regression, taking into account that 
each star ; contributes to v\\ with its own error Avy,- (referred to 
as Av^ AVE (/,-, bi) in Eqn. (|6]l). Hence, for instance, the error on 

the mean will be Av 2 = (P + ) Av 2 , on the dispersion tensor it 

2 

w ill be Act 2 = 2((P <g>P) + ) (v'y ® Av',,) 2 , and so on (exactly as 
in lFuchs et ail J2009I) '). 




GC 



Fig. 2. Schematic representation of the hypothesis of an 
anisotropic and isothermal distribution of an unbiased sample 
of thick disk velocities superimposed on a thin disk sample. The 
observer is located at the Sun's position R G at the centre. The 
Galactic centre (GC) is at the bottom of the figure. The direc- 
tions on the celestial sphere to each star s, are indicated with 
dashed arrows for thick disk stars and solid arrows for the thin 
disk stars. The green thick arrows represent the radial velocities 
of the thick disk, whose corresponding unknown true velocity 
vector is in light blue. The thin disk radial velocities are repre- 
sented by the thin green arrows, and the corresponding unknown 
olive-green arrows represent the true velocity vectors. Note that 
the olive-green arrows are generally longer in the GC direction 
and shorter in the anti-Galactic centre direction as expected from 
the velocity dispersion trend. 



3.2. The scattering processes and the thick disk component 
disentanglement 

So far, we have a procedure that, through the Eqns. ( ITOl . 
( fTTT i and ( TT2b . provides us with the first four true moments (i.e., 
fy, er, S, T}) of the distribution function of the thin and thick 
disk mixture from the projected moments along the l.o.s. (i.e., 
{vy, <T||, Sy, T||}) that we obtain directly from the observations. 

To proceed further with the thick disk analysis, we want to 
disentangle the two distribution functions (DFs) and proceed by 
determining the thick disk velocity dispersion tensor alone. 

We assume that the two single-particle DFs for the thin and 
thick disk, fa (x, v) = / t hi n and fa (x, v) ee / t hi c k, can be added lin- 
early to yield the overall distribution function / tot (x, v) = fa +fa. 
For the mixture distribution function / tot in the previous Section 
we computed the first four moments. The first and the second 
moments correspond to the mean of the sample and to the dis- 
persion tensor, respectively, 

cr 2 = ±j Jy'^f M (x,y)dy, (13) 
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where the integral extends over the whole velocity space and N 
is the number of stars. The third and fourth moments are: 



and 



= jff V 8 7 ,ot (x,vMv, 
= v' 84 / tot (x,v)rfv, 



(14) 



(15) 



with the same definition interval for the integrals. We see that 
N - N(x,t), as well as er 2 = er 2 (x, f), S = S (x, t) and 
T = T (x, t) because, in general, f tot = / tot (x, v, f)- The form 
of this distribution function is generally unknown but its sub- 
components, / t hin or /thick, are often assumed to be quadratic 
in the peculiar velocity components, i.e., of a generalized 
Schwarzschild type. 

In the time evolution of the total DF / tot , we can argue 
that a very different role is played by /hi n or /hick- A few pro- 
cesses of scattering are generally assumed to influence the time 
evolution of /hi n : mostly the scattering due to transient spi- 
ral arms and the encoun ters with giant molecular clouds (e.g., 
Binn ev & Tremainell987l) . Both processes lead to an increase of 
the velocity dispersion with increasing distance from the plane 
and provide insight into the origin of the age-velocity dispersion 
relation. Moreover, the Schwarzschild characterization of / t hi n , 
on large scales in the velocit y space, is not as good a d escrip- 
tion on smaller scales (e.g., Seabroke & Gilmorel l2007l) . Sub- 
structures are mainly found for young stars (but not only for 
these). Such stars in a moving group are born at the same place 
and time, and then disperse into a stream that may intersect the 
solar neighbourhood. In this scenario, the stars that are moving 
in the same group shoul d share the same age, metallicity, and 
azimuthal velocity (e.g., iDehneri [19981 IChereul et ai]ll998l) . A 
different e xplanation for the sub-s tructures in the velocity DF 
is given by De Sim one et al.l (120041) . These authors suggest that 
sub-structures arise naturally from the same spiral gravitational 
fluctuations that excite the growth of the velocity dispersion. In 
this picture, sub-structures are caused by homogeneous star for- 
mation in an irregular potential, as opposed to inhomogeneous 
star formation in a regular potential. Once the thin disk compo- 
nent is disentangled from the thick disk component, its analy- 
sis will proceed with the full exploitation of the proper motions 
(Paper II). 

The situation is less clear for the thick disk where the role 
of the scattering processes is still undetermined. Its formation 
may be related to the influence of massive satellites that can 
either heat pre-existing disks (preservin g a vertical metallic- 
ity gradient if present. (|Katz etal.N2011l) or get accreted (e.g., 
lOuinn & Goodman! Il986t lOuinn et al.lll993h. Other formation 
scenarios include gas-rich mergers (e.gTlBrook et al. 2004), or 



intern al radial migration processes (e.g., Schonrich & Binnev 



2009) or formation induced by satellites (e.g., Di Matte o et al 



2011). We will start working with the hypothesis that the thick 
disk DF can be approximated as /tuck K e ^ x ' c " ,ilkV ') where 



the covariance matrix 



cr 1 

'-'thick 



is er thick and the quadratic equa- 



tion v' r C t hickv' = 1 defines a velocity ellipsoid that character- 
izes the distribution function of the thick disk fully and that we 
want to determine. We will check "a posteriori" (see Section 
14.4.2b if this approximation is consistent with our results for 
the group of stars selected as in Section 12.31 The thick disk is 
more diffuse in the solar neighbourhood tha n the thin disk and is 
chara cterized by a larger scale height (e.g., Binn ev & Merrifieldl 
1998). We assume that the dependence of the velocity ellipsoid 



in the sample analysed with RAVE can be adequately consid- 
ered to be independent from the Galactic plane distance, i.e., an 
isothermal but anisotropic DF, as commonly assumed in most 
of the thick disk velocity ellipsoid determinations available in 
the literature (e gjLayden et alJI 19961: 1 Alcobe & Cubarsill2005l: 
IChiba & Beersll2000t iHolmberg et al.ll2007l) . We will return to 
the validation of this hypothesis in Section l4~4l 

Under these assumptions, the decoupling between the thin 
and thick disk co mponent is straightforward and analytical ly de- 
veloped in, e.g. JCubarsll dl992h (see also ICubarsil d2007l) V We 
are going to ap ply this inversion techniqu e in its version based on 
the cumulants ( Cubarsi & Alc obe (2004)) because this is numer- 
ically more stable and hence suitable to minimize the numerical 
error propagation. 

We outline the methodology applied here because it differs 
slightly from the original work o flCubarsi & Alco be ( 2004). The 
interested re ader is referred to Append ix B, where the technique 
contained in ICubarsi & Alcobei d2004l) is explicitly derived and 
adapted for our specific case. The derivation of the cumulants 
from the Eqns. d!3l . (TT4b and ( fT31 ) is straightforward: 

N 



N- 1 



N 1 



(N-l)(N-2) 

N 2 (N + 1) 
(N - l)(N - 2)(N - 3) 



(16) 



N+l 



hereafter simply indicated as Kijk, kiju, respectively. If we let 
Nam and Mhick be the unknown numbers of stars for the thin and 
thick disk components, we can define the parameter q as follows: 



q = 



»thick 
"thin 



n thin 
"thick 



(17) 



where « t hi n = ^ thin /^f and « t hick = ^ thicl ^v are tne unknown frac- 
tions of stars in the thin and thick disk. This parameter varies 
from when 50% of the stars belong to each component to oo 
when 100% of the stars belong to the thick disk. Our goal is 
to maximize this parameter q in the subsample of stars selected 
fro m the RAVE catalogue as in Section [2~3l 

ICubarsi & Alcobei d2004l) show the great advantage of ex- 
ploiting the symmetries of the normal distribution functions by 
working with a transformed peculiar vector. We follow this con- 
vention by defining a normalized vector as 



d = V n thin«thick > 

W2 



(18) 



where w = {w\, W2, W3] is the velocity vector of the difference 
between the thin and thick disk centroids; w = v t hi n - v t hick- The 
reason for the normalization to the highest velocity difference 
component of w is to reduce the error propagation in the SVD 
solution using the best-determined differential velocity, i.e., W2- 
We define the transformed peculiar velocity V that keeps the dif- 
ferential velocity in the azimuthal direction invariant: 

V = H v 



H 



-d 3 r 
1 

-1 d x 0, 



This is an isomorphic transformation of the peculiar velocity 
vector v, with det (H) = 1 due to the previously assumed normal- 
ization of Eqn. dT8) , The new cumulant equation for the trans- 
formed DF function of V can be derived as a function of the 
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old cumulants of Eqn. ( fT6b (see Appendix B). Once the new 
V-mixture cumulants have been computed, say K [3] and K [4] 
for the third and fourth cumulants respectively, we determine 
the parameter q of Eqn. (flTT i by making again use of the SVD 
technique to determine the solution of the linear system of over- 
determined equations as already done before for the system of 
equations in Section |3~T1 We present the description of the sys- 
tem in Appendix B to which we refer the interested reader. 

Finally we extrapolate the values of «thick, ""thick, and v/%, 
i.e. the characterizing thick disk kinematics parameters we were 
searching for, with a simple numerical iterative procedure that 
maximizes the thic k disk componen t alone (e.g., Simulated 
Annealing Methods, iPress et al.l 1 19861) . Within the preselected 
set of data satisfying the selection criteria of Section l2~3l we con- 
sider a Monte Carlo generation of the stellar parameters within 
their observational errors and we isolate the thick disk compo- 
nent alone with the methodology described above. The kine- 
matic characterization of the thick disk is achieved with the same 
purely geometrical approach developed for the mixture but now 
on the thick disk selected sample alone. Once we apply this 
Monte Carlo approach to generate several catalogues, and then 
apply the procedure to each of them, we summarize our results 
after a simple statistical analysis. 

4. Results 

Here we present the results for the determination of the motion 
of the Sun relative to the LSR and the velocity dispersion tensor 
for the thick disk. 



4.1. Solar motion relative to the LSR and thick disk rotation 
velocity 

Computing the first moment of the velocity part of the distribu- 
tion function, namely the mean, we are left with a free param- 
eter ||v f |. We assume that the dependence of ||v c || on the model 
is rather weak for the two components, v^o and Vw,®, which we 
hence derive as 



Vu,e 

Vf,0 



(9.87 
: (8.01 



: 0.37) km s" 
: 0.29) km s" 



(19) 



The determination of the last value of v Q : vy, is more compli- 
cated. This is mostly for the following two reasons: 

1 . The overlap of the rotational delay of the population vi ag and 
the contribution of the peculiar motion of the Sun in the same 
direction Vv,e, both acting together on the mean streaming 
velocity ||v c ||, prevents us from obtaining the relative veloc- 
ity between the LSR and the standard of rest centred on the 
motion of the Sun from our first distribution moment, \. Thus 

remains undetermined. 

2. A full determination of \ Q could in principle be achieved 
with this methodology by including all stars down to a low 
latitude, say, b < 15 deg. However, the survey data are miss- 
ing this low-latitude sky coverage. Moreover, there is no 
mapping in the Galactic rotation direction (see FigureQ}. 

Nevertheless, while we cannot carry out a clean determina- 
tion of the first moment of the velocity section of the DF, by 
assuming the literature va lue vy i0 = 13.5 ± 0.3 km s from 
iFrancis & Anderson! d2009l) . we gain a value for the velocity lag 
of our thick disk sample of stars, 



Table 1. Values of the thick disk velocity dispersion tensor. 
Absolute values of the square root of the off diagonal elements 
are reported. 





[kms" 1 ] [kms" 1 ] [kms -1 ] 


°~ RR °~ R<t °~ Rz 

a zz 


59.2 ±4.4 31.1 ±21.2 9.3 ± 2.3 
47.3 ±7.5 3.2 ±7.1 
35.9 ±4.1 



As explained at the end of the previous section, the procedure 
works iteratively in order to select a single population of thick 
disk stars from which to deduce the kinematic parameters. This 
is done because we cannot determine the solar motion relative 
to the Solar LSR from a non-local mixture because the thin disk 
velocity dispersion trend on the meridional plane could bias the 
results. Thus we need to disentangle one single component from 
which to derive (vu,e, vw.e)- 

4.2. Thick disk velocity ellipsoid 

The values determined for the velocity dispersion tensor, ctrr = 
59.2 ± 4.4 kms -1 , cr 00 = 47.3 + 7.5 kms -1 and cr zz = 
35.9 + 4.1 km s (see Table Q]), are in good agreement with 
what was alread y presented in the literature for the thick disk 
component (e.g., Alcobe & Cubarsi 2005; Chiba & Beers 2000; 



iHolmberg et al.ll2007 ). The presence of the thin disk is expected 
to bias these values only slightly in the sense that the random 
selection of the stars used in order to maximize the parameter q 
in Eqn.dTTll is proven to be able to sep arate the two components 
as shown in ICubarsi & Alcobd d2004t) and as confirmed by the 
simulated mock catalogues (see Section l4~3l ). 

As can be seen, the vertical tilt of the velocity ellipsoid is de- 
duced from the mixed term component cr Rz = 10.1 ± 3.3 kms 1 , 

leading to a tilt angle of about = j arctan ^ ^J^ r ) — 3 -07 ± 

1.1 deg. This is similar to the result of Siebert et al. (2008), who 
applied different selection criteria to the same data, but deviates 
from the recent determination bv lFuchs et al.1 (120091) . Regarding 
the comparison with these studies it is worth noting that we used 
Eqn. d6}, while in lFuchs et al.l d2009h the tilt depends on the mix- 
ture of the two components (thin and thick disk) and it was con- 
sidered without any technique capable of disentangling the rela- 
tive statistical influence of the two, nor accounting for the radial 
dependence on R in cry = o-jj(R,z). In the meridional plane the 
thin disk's radial component profile could partially influence the 
determination of the velocity ellipsoid tilt, if the binning is done 
only in the vertical direction (as proposed in their paper) but not 
in the radial direction. 



4.3. Testing the results with the Padua Galaxy model 



Vthick,i aK = 49 ± 6 km s 



(20) 



As indicated at the end of Section 13.21 in the spirit of a Monte 
Carlo approach to the error analysis, we proceed by generating a 
set of synthetic catalogues, numbered from j = 1, iV ca t (with 
A^ cat being a high number) by assigning to the i th star a ran- 
dom distance d, within its photometric distance error Ad, and 
a random radial velocity v n ; within its radial velocity error Av r ,, 
(no errors are assumed in the star's coordinates (I, b)). Then, the 
procedure to disentangle the stellar thin and thick population is 
applied for the / h realization of the catalogue and for all the 
j = 1 , A^at. If A^ C at is sufficiently high, the stability of the error 
is achieved with standard statistics tools on the N C3l realizations 
of the catalogue. 
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Nonetheless, to be sure our method is working correctly 
and to get more insight on its limits and potential, we apply 
the methodology just described to a fully analy tical catalogue . 
Then we apply the P adua Galax y model (e.g. , iNg et al.l l 1995. 
I2002t I Vallenari et all 120061) and iPasettol (l2005l) and references 
therein) to generate a full synthetic Galaxy where everything is 
under control to test the ability to recover the properties of the 
data generated by the new method developed here. In this pro- 
cedure we proceed by producing mock RAVE catalogues with 
prescribed relations and parameters that we ultimately want to 
recover with the technique here developed. Initially three syn- 
thetic stellar populations, representative of the thin and thick 
disks and halo are generated. T he s tellar structure mode ls are 
taken from iBertelli et all d2009l) and iBertelli et all (120081) from 
which we generate a stell ar population by as suming an initial 
mass function (IMF) (e.g.. lKroupa et al.ll 1993h and the star for- 
mation rate (SFR) is assumed as a free parameter. Once the syn- 
thetic Hertzsprung-Russel diagrams with N stars of the stellar 
populations are generated and projected in the corresponding 
RAVE/2MASS colours and magnitude passbands, we can dis- 
tribute the N synthetic stars representing the mixture of three 
stellar populations in the phase space of the M W once a model 
of int erstellar extinction is assumed (e.g., from Drimme l et alj 
120031 We are particularly interested in the kinematic descrip- 
tion the stellar populations of the disks. For them we assume 
double exponential spatial density profiles 

where t is age of the stellar populations that we discretize for 
simplicity in four temporal ranges for the thin disk: i = 1 for t e 
[0, 3[ Gyr with p ,t = 2.2 x 10 8 M G kpc~ 3 , h RA = 2.90 kpc, h z ,i = 
0.20 kpc, i = 2 for t e [3, 5[ Gyr with p , 2 = 1 .0 x 10 8 M G kpc~ 3 , 
h R>2 = 2.90 kpc, h zA = 0.25 kpc, i = 3 for t e [5,7[Gyr with 
p ,i = 2.1 x 10 8 M G kpc~ 3 , h Rtl = 3.10 kpc, = 0.28 kpc, i = 4 
for t € [7, 10[Gyr with p ,i'= 8.0 x lO 8 M kpc~ 3 , h RA = 3.10 
kpc, h Zi \ = 0.35 kpc. The thick disk scale parameters adopted are 
Po,Mck = 1-1 x lO 6 M kpc~ 3 , h R ^i = 2.50 kpc, h z ^ - 1.31 kpc and 
interstellar medium parameters are po,i = 1.48 x lO 8 M kpc~ 3 , 
hu.i = 4.54 kpc, h Zt \ = 0.20 kpc. To complete the density profile 
from which we deduce the potential that is used in the Boltzman 
equation to implement the kinematic description (as, e.g., in Eqn. 
(l24l >) we add a simple Hernquist bulge (Hernquist 1993) and log- 
arithmic potential for the halo. We then tune a Poisson-solver to 
match the observational constraints on the rotation curve, i.e. the 
Oort functions, the local density ratio between the stellar popu- 
lations, the terminal velocity for the inner Galaxy with respect 
to the solar radius, the tota l mass inside 100 kpc etc. (see ap- 
pendix A in Valle nari et alj d2006l) for an extended description 
and Pasetto 2005). 

It is beyond the goal of this paper to derive the structural 
parameters of the MW from the RAVE data. Thus, the previ- 
ous values are frozen and have to be considered as guess values 
assumed to test our new methodology. Other literature scale pa- 
rameters, short er in the scale length of the thick disk, c an be eas- 
ily tested (e.g., ICarolloet ai1l2010t iBensbv et al1l201 lh as sug- 
gested in relation to the alpha-enhancements (e.g.. lCheng et al.l 
2012; Brook et al.ll2012l) . Once we have fixed these parameters, 
we fix the trend of the velocity dispersion tensor for the thin 
disk in the meridional plane crj^ (see Paper II) and arbitrary 
values for the velocity dispersion tensor of the thick disk 0"jL ck 
with which to generate our mock catalogue and that we want 
to recover with our novel method. The Padua Galaxy modelling 



technique produces a mock catalogue with N ^ N synthetic val- 
ues for 

{m,c,logg,r hel ,l,b, v r , p,} (22) 

where m (and c) are magnitude (and colour) in a desired pass- 
band, and \og l0 g is the logarithm of the stellar surface grav- 
ity. Considering that the cut adopted in the surface gravity in 
Section 12.31 lies around the turn-off of a thick disk population 
with t e [10, 12[Gyr, we have an almost bijecti ve relation be- 
tween stellar magnitude and distance once the iDrimmel et al.l 
(2003) 3D extinction map is considered, thus rhei is a catalogue 
of photometrically determined distances for the dwarf stars of 
the RAVE catalogueQ In Eqn. d22]l, when N - N, then /, b are 
exactly the observed RAVE stellar directions, v r (and /}) are the 
radial velocities (and eventually proper motions) for our synthet- 
ically generated catalogue. 

We extensively tested the new technique against the ability 
to recover 

- the correct thick disk velocity dispersion tensor with which 
we generate the RAVE mock catalogues, 

- the influence of large-scale effects on the thick disk velocity 
ellipsoid (see also Paper II), 

- the minimal sky coverage for which we can correctly recover 
the underlying thick disk kinematics, 

- the possibility to correctly disentangle the Vy iQ > 

- the role of the log 10 g to disentangle dwarf and giant stars in 
our sample. 

But one of the most interesting results of our extended analysis 
on the stability and selection cut criteria is probably the exercise 
presented in Appendix C, where we show that the methodology 
developed here mathematically can retrieve the correct thick disk 
velocity ellipsoid even without the use of photometric distances. 

It is possible to make further cuts, for ex ample in the colour 
magn itude diagr a m (CM P), but as shown in I Veltz et al.l (|2008) 
or ISiebert et alj (|2008) this drastically reduces the available 
number of stars. Large numbers of stars, however, are the means 
of our present method to reduce the error bars. Interesting dif- 
fer ent approaches hav e also been discussed i n t he literature, e.g. 
bvlRobin et alJ (l2003l):ISchuster etall (Il993l) or lCasertano et al.l 
( 1990); Ratnatunga et al. ( 1989|). The applicability of their meth- 
ods to RAVE data is also interesting but beyond the scope of the 
present paper. 

4.4. Discussion of the results 

In the following subsections we will first take a look at the as- 
sumptions and approximations trying to verify them aposteriori 
and then we will improve our results. 

4.4.1 . Looking to the past 

Up to now, we approached the problem of the determination of 
the thick disk mean velocity and its velocity dispersion tensor. 
This approach permitted us to determine the kinematic proper- 
ties of the thick disk component in the framework of an isother- 
mal velocity ellipsoid model (Table[TJ. Nevertheless, in the liter- 
ature the ellipsoid is commonly described as a velocity distribu- 
tion independent of Galactic position because of the difficulties 
in achieving a more accurate description. However, both these 
hypotheses are without a rigorous theoretical basis for the thick 

1 Apart for the test performed in this section, we will adop t every - 
where the photometric distances as determined in lZ witter et all (2010). 
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Fig. 3. Plot of the errors in velocity and position: Avf AVE is 
the error in the radial velocity as observed in the RAVE sur- 
vey. Arhei is the erro r in the photometric distances as inferred in 
Zwit ter et all d2010t> . The binning is 0.5 km s 1 and 50 pc with 
a scale colour normalized to the highest number of stars per bin 
(in light yellow) to stars per bin (dark blue contour). 



disk component of our Galaxy. Here, thanks to the RAVE survey, 
we can investigate this aspect further. 

In Figure [3] we plot the errors in the velocity and position for 
the sample of stars in the mixture of thin and thick disks as ob- 
tained from the catalogue (same stars as in Figure[T). As already 
pointed out in Section 12.11 we work with the stars with larger 
errors in the velocity. Thus the bulk of the sample has radial ve- 
locity uncertainties up to 6 km s _1 , slightly hig her than the av- 
erage error produced by the RAVE survey (e.g., Steinmetz et alJ 
2006). The photometric distance err ors are typically of t he order 
of 0.15 kpc for our selected sample (Zwitter et alJl2010t) . 

A free parameter of our modelling approach is the mean az- 
imuthal velocity of the RAVE sample we select. This quantity is 
related to the asymmetric drift of the stellar populations, the ro- 
tational lag of the thick disk and the circular velocity of rotation. 
The thick disk subsample of the RAVE catalogue analysed here 
appears to lag behind the Solar LSR by = 49 km s . The gen- 
eral behaviour of the asymmetric drift is not known a priori, es- 
pecially in its vertical trend and for a mixed set of thin and thick 
disk populations. For a single stellar population in the disk, the 
mean circular streaming velocity can theoretically be predicted, 
starting from the Jeans equation in cylindrical symmetry. We can 
approximate the term \ c (x) in Eqn. (Q]) as 



v c (x) = \ c (R, z) 



( \y c (R,z)\ r ™ c ™ bsin > > 

\v c (R,z)\ R °- rb ° > ™ bcosl 
10 



with 



ivc(^ Z )i = K S r-I^Kr+0 +( 

T n ^ dR dz dR j 



n 4r + 

1/2 



(23) 



(24) 



Table 2. Overview of the data sets and type of analysis per- 
formed by recent literature works (not based on RAVE data). 
Column one refers to the number of stars effectively used in 
the analysis of the thick/thin disk, not to the total number of 
stars produced by the database. HR/LR in column two stands for 
high/low spectral resolution respectively as specified in the pa- 
pers of column three. For comparison we report the number of 
stars used in the present study in the last row. 



no. stars 


Type of analysis 


References 


451 


Chemistry (HR spectra), kinematics 


(1) 


176 


Chemistry (HR spectra), kinematics 


(2) 


250 


Chemistry (HR spectra), kinematics 


(3) 


1,498 


Strmgren photometry, kinematics 


(4) 


306 


Chemistry (HR spectra), kinematics 


(5)(2)(9)(11) 


412 


Kinematics (HR spectra) 


(6)(7)(12) 


76 


Chemistry (HR spectra), kinematics 


(8) 


102 


Chemistry (HR spectra), kinematics 


(9) 


17,277 


Chemistry (LR spectra), kinematics 


(10) 


23,767 


Chemistry (LR spectra), kinematics 


(13)(14) 


38,805 


Kinematics (HR spectra) 


present work 



Refe rence s 

Fuhrmann 



(1) lNevesetal.1 ([2009]); (2 ) iReddv et all (120061); (3) 

J2004h: (4) ICasagrande et all 120111): (5) iNavarro et al.1 

j2011h:(6)[M oni Bidi n et al.l d2012tJ); (7)lMoni Bidin et al.||2012bl) : (8) 
Afsar et al.ld2012h: (9)lBensbv et all (120051): (IQ)lLee et al.U2011|): (11) 
Nissen & Schuster! d2010h: (12) iGirard et all ( 120061) : (13) iBovv et al l 
<2012ah : (T4) lBovv et al . [2012b). 



and the density profile act together to play a crucial role in the 
prediction of the streaming circular velocity and in its axisym- 
metric approximation. Assuming that the selected data sample is 
mostly representative of the thick disk component, the previous 
Eqn. d24b is reduced, to: 



( d\np(R ,O) 
i din/? 



+ 1 



RR 



HL + vl = 201 km s" 1 



where Oxot is the total potential of the MW. This indicates how 
the terms of the velocity dispersion profiles, the circular velocity 



where /i« = 2.9kpc is the adopted value for the scale length o f 
a double exponential density profile, e.g.. lVallenari et all 0004). 
v c is the circular velocity of the Galaxy in the solar neighbour- 
hood as above, and <j 2 RR and <x^ take the values derived in our 
paper. The expected value for a pure thick disk component is 
not so far from the best fit value derived in our study, namely 
||v c || = 179 km s . We point out here that this permits us to im- 
pose only an upper limit on the expected value of the streaming 
circular velocity for the mixed sample of thin and thick disks we 
analysed. This is thus not necessarily in dis agreement with pre- 
vious studies predicting lower values (e.g., Soub iran et al.f 2003) 
because of the thin disk influence. A further cut in the latitude 
can in principle reconcile these values but unfortunately reduces 
the number of stars that we can retain with our iterative pro- 
cedure and consequently increases greatly the error bars in the 
velocity dispersion tensor. 

Thus we prefer to consider ||v e || as a free tuning parameter in 
our approach, remembering that if the thick disk results from 
a sudden heating due to an infalling satellite, it can be kine- 
matic ally more decoupled from other components of the MW 
(e.g. lOuinn &~ Goodman Tl986tlWalker et alJll996l: ISpagna et al.1 
2010; Bekki & Tsuiimoto 201 1). The overlap between rotational 
velocity and peculiar motion of the Sun relative to the LSR can 
not be disentangled with the use of an isothermal, non-local sam- 
ple of sta rs so it remains beyon d the scope of the present work 
(see e.g.. lSchonrich et aD2010l) . 
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4.4.2. Looking to the future 

In principle, it is nowadays possible to disentangle the thick and 
thin disks using selection criteria based on the chemical pr oper- 
ties o f the stars (e . g., see the recent review by lNissenll201 ll) . In 
brief, iNeves et al ] d2009h determined the trend of abundance ra- 
tios as a function of [Fe/H] from 451 high-resolution spectra of 
F, G, and K main-sequence stars in the solar neighbourhood con- 
firming the long known bimodal distribution of the [a/Fe] ratio 
for the disk stars (with the thin disk stars less alpha- enhanc ed 
than the th i ck dis k ones), see for instance Ikeddv et ail d2006l) or 
iFuhrmannl (120041) and references therein. Along the same line of 
thought. [Casagrande et al.l (1201 ll) determined the ratios [Fe/H] 
and [a/Fe] for 1498 selected stars of the Geneva-Copenhagen 
Survey from Stromgren photometry, even if their data do not 
provide a clear bimodal distribution betwe en thin and thick disk 
stars as far as the [a/Fe] ratio is concerned. lNavarro et al.l d201 lb 
investigated the separation of thin and thick disks using a com- 
bined index of [Fe /Hl, [g/Fel and the heavy element Eu for a 
sample of 306 stars. iRuchti et al.l (1201 lal) studied the vertical and 
radial gradients in metallicit y and alpha-elements for a sample 
of selected thick disk stars. iMoni Bidin et al ] (I2012allbh studied 
412 red giant stars in the direction towards the South Galactic 
Pole trying to describe kinematics, chemist ry, and content of 
dark matter in the MW di s ks (bu t see e.g., iBovv & Trem aine 
120121) . Finally, lAfsar et al.l (1201 2|) conducting a detailed abun- 
dance analysis and atmospheric parameters of 76 stars in the thin 
and thick disks with very high-resolution spectra (^R ^ 60, 000) 
reaching conclusions similar to those found by other authors. 

Common to all the above studies based on chemical prop- 
erties (see Table 14.4. U . is the small number of stars that are 
considered and the limited spatial coverage in about the solar 
neighbourhood. In contrast, our method, owing to its statistical 
nature, first requires a high number of stars such as that produced 
by RAVE, and second does not need chemical information to 
work properly. It is worth recalling that data on radial velocities 
(from medium resolution spectroscopy from which the kinemat- 
ics is derived) and photometry are easier to gather than good 
chemical abundances (from high resolution spectroscopy). Our 
approach can easily be applied to large sets of spectroscopic data 
such as the continuously growing RAVE survey (or Gaia in the 
near future). Finally, our method makes use only of kinematic 
data, leaving the coupling between kin ematics and chemis try to 
be investigated in a forthcoming paper dBoeche et al.ll2012h 

Although in recent years we have seen continuous im- 
provements inGalactic modelling, see e.g. the Besanon model 
(Robin et al. 2003), the Pado va Galaxy m odel (used in this pa- 
per), DF based mod els (e.g.. lBinnevl20 1 2h or the Galaxia model 
dSharma et al.ll201 ll) . we still lack a good theoretical framework 
coupling kinematic s and dynamics with population synthesis 
and chemistry (e.g., IScho nrich & Binne^HOQa I Just et alJl20lU 
Uust & Jahreil2010l) . 

Neglecting the radial and vertical dependence for the thick 
disk velocity ellipsoid is an almost universal assumption that in- 
deed produces acceptable results for the vertical tilt in our case 
as well. Nevertheless, we can further explore this working hy- 
pothesis by making use of the available photometric distances. 
Although the proximity of the stars selected grants the validity of 
our approximations, it limits the exploration of the more distant 
zones of the thick disk. For every sample statistically represen- 
tative of the thick disk population that we obtain with the pre- 
viously outlined method, we split the data set into two regions: 
within and beyond the Sun's Galactocentric position, R < R Q and 
R > R Q respectively (where R is the cylindrical Galactocentric 



Table 3. Values of the thick disk velocity dispersion tensor for 
R < R Q . Units as in TableQ] 





[kms" 1 ] [kms" 1 ] [kms~'] 


°~RR 0~R(/, C Rz 
£7" 0j, 
(T- z 


60.2 ±7.1 37.6 ±21.7 13.3 ±9.8 
44.7 ±8.1 4.0 ±7.2 
37.2 ±5.7 



Table 4. Values of the thick disk velocity dispersion tensor for 
R > R Q . 





[kms" 1 ] [kms" 1 ] [kms" 1 ] 


CRR 0~Rp a Rz 

tr <i,z 

(T-, 


55.8 ± 6.5 35.5 ±20.2 9.6 ±7.6 
45.2 ±7.3 3.8 ±3.1 
36.3 ±4.1 



radius). The values obtained for the velocity dispersion tensor 
are then statistically averaged, as done to obtain Table[TJ and are 
listed in Tables [3] and [4] 

The entries of Tables[3]and|4]hardly show any difference be- 
tween the stars outside and inside the solar circle. Nonetheless, 
the suspicion arises that within the errors, there may be some 
evidence for a general increase in the velocity dispersion of the 
thick disk along the radial direction (the velocity dispersion is 
higher for star inside the solar circle than outside). This implies 
that for the thick disk alone the classical isothermal picture needs 
to be improved. The thick disk velocity dispersion tensor seems 
to imply a dependence on the position in the meridional plane 
^thick _ p-thick^j z ) rj-jjg e ff eC ( j s> however, small due to the 

small range of distances sampled with our data. We have higher 
error bars in the inner sample because of the cut we have applied 
on the radial velocity errors (see Eqn. (O) that retains a higher 
number of stars in the anti-centre direction (see Figure [TJ. As 
a consequence, despite the stellar density decrease towards the 
outer regions of the Galaxy, the error bars of the inner and outer 
samples of stars are comparable in size. 

This suggests a way of evaluating the effect of the dis- 
tance if a gradient in the vertical dependence of the absolute 
value of the azi muthal velocity occurs along the thick disk (e.g., 
Chib a & Beersl2000T) . Eqn. (|24"jl. which hold s for small distances 
from the plane (e.g.. IVallenari et all 12006, their Eqns. (31) to 
(33)) can be used to create mock catalogues with radial and ver- 
tical gradients in \\' c (R, z)\ where i refers to the thin and thick 
disk once the populations are embedded in a common potential 
O tot . In this way we can check the role of small gradients, if any, 
in cr^ ick = of? ck (R,z) on |v* ick (R,z)| or we can simply artifi- 
cially place a gradient in |v* lck (/?,z)|. This approach has been 
tested with Eqn. (l24l i and gradients up to 30 km s _1 kpc _1 . The 
resulting values for the thick disk velocity ellip soid are withi n 
the error bars presented in Table[T](see also Be rtelli et al.l l2003). 
On the one hand, this is because of the limited vertical and ra- 
dial extension reached by our sample, which reduces large scale 
effects. On the other hand, this happens due to the method itself 
not involving a selection based on chemical abundances. For ex- 
ample, a vertical gradient in mean azimuthal velo city of about 
20 km s kpc was found bv lSpagna et al.l (1201 Oh even though 
their result refers to \z\ 6 ]l,3[kpc, i.e., a region beyond the 
range spanned by our sub-sample of th e RAVE c atalogue. 

The different result in the work o f IVeltz et al.l J2008) is due 
to their data sample being based on fewer than 600 stars in 
the direction of the south Galactic pole (SGP). This leads to a 
higher value of cr RR , on which the influence of halo stars is yet 
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to be investigated (see, e.g., IChiba & Beers 2000; Smith et aT] 
2009). Particular attention should be paid to the interpretation of 
the components that are mixed with the azimuthal component: 
cr S0 = 29.4 ± 17.2 km s~\ cr 0z = 5.8 ± 5.1 km s -1 . The first is 
reminiscent of the well-known vertex deviation. For the thin disk 
the vertex deviation can be predicted from ctrq and amounts to 
an angle of = 21°. Its large error is probably due to the thick disk 
components closer to the MW plane which could induce circular 
velocity gradients that we were unable to detect with our method. 
The <70 Z component is a further indication of the coupling of the 
vertical and azimuthal velocity components. 

Finally we point out that the method we have developed can 
be extended to an arbitrary number of populations provided that 
first the formalism is expanded to include higher order cumu- 
lants, second each population is suitably sampled (sufficiently 
high number of stars), and third each population is characterized 
by at least one distinct kinematical parameter. In our case, we 
present evidence of the existence of at least two distinct popula- 
tions whose characterizing parameters are the second order ve- 
locity dispersion tensors. To prove this statement we perform the 
following experiment. Suppose that the solar position is moved 
closer towards the Galactic centre, for instance at the position 
Rq < 5 kpc. Then we set up a mock catalogue centred on the 
new position of the Sun, R Q . By construction, the velocity dis- 
persions of the thin and thick disk stars are the same (no cor- 
rection for different extinction is applied). Since the kinematical 
parameters of the populations are identical, the method fails to 
converge to a solution. This result can be better understood when 
looking at Fig|2l where now we imagine that the Sun is located 
at R Q and the vectors representing the velocity distribution of 
thick disk stars (light blue arrows) have the same length of those 
for the thin disk stars (olive-green arrows). Therefore, there is no 
kinematic way of distinguishing two different populations from 
the radial velocity distribution of the mixture (the thick green ar- 
rows). However, this does imply that the two populations could 
be separated by considering other parameters such as chemistry, 
alpha-enhancements etc. 



5. Discussion and conclusions 

In the sample extracted from the RAVE catalogue that we anal- 
ysed, two populations of stars with different kinematics are 
found that correspond to the thin and thick disks, even though 
our analysis actually focused only on the thick disk. 

The major difference with respect to previous studies in the 
literature is the peculiar spatial volume covered by RAVE, which 
allows us to determine the kinematics of the thick disk not by 
extrapolating from a single small field-of-view to the whole thick 
disk, but by directly measuring the kinematics by studying the 
much larger volume of the thick disk observed by RAVE (see 
also Paper II). Over this extension, the key thick disk parameters 
are determined as: 

1. two components of the solar motion relative to the solar 
LSR, namely as v^ Q = (9.87 ± 0.37) km s and vw,o = 
(8.01 + 0.29) km s -1 ', 

2. the rotational lag of the thick disk component relative to the 
LSR v t hick,ia g = 49 + 6 km s" 1 , 

3. the velocity dispersion tensor of the thick disk considered 
to be an isothermal population: cr RR = (56.1 ± 3.8)kms _1 , 
cr Rip = (29.4 ± 17.2)kms _1 , cr Rz = (10.1 ± 3.3)km s -1 , cr w = 



Moreover, we mention that the missing full determination of 
the Sun's velocity vector relative to the LSR is just a choice. 
This does not mean at all the impossibility for RAVE to charac- 
terize completely the solar nei ghbourhood (see, e.g.,|Kiss et al. 
120101: ICoskunoglu et al.ll2011l) . We defer the study of the thin 
disk kinematics to a companion paper (Paper II). 

The determination of the thick disk velocity dispersion ten- 
sor and its behaviour in the (0,R,z) plane is just a small step 
in the investigation of this Galactic component. The presence 
of a small gradient, especially in the cr RR component, does not 
represent the failure of the classical picture of an isothermal de- 
scription within the range of distances investigated but it is an 
example of the quality of data that a radial velocity survey like 
RAVE can provide to confirm, extend or investigate new ideas. 

The thick d is k is a promin e nt fea tu re of our galaxy 
(iGilmore & Reidl 119831: lYoshiil 119821: I.Turic et al.1 120081: 
l Yannv et al.l 12009}) and of e xternal disk galaxies (e.; 
van der Kruit & Fr eemanI 1201 It lYoachim & Dalcantonl [2006 
Several possible formation mechanisms have been suggested 
for the MW thick disk formation. For instance it may be 
related to the influence of massive satellites that can either 
heat pre-existing d isks or contribute by being accreted (e.g ., 
l Ouinn & Goodman[ll986tlOuinn et al.lll993tlAbadi et al.l l2003). 
Wal ker et al.l 11996) showed in detail how low-mass satellites, 
while rapidly sinking into th e potential w e ll of a galaxy, could 
substantially heat a disk. Br ook et alj d2004l) investigated 
the influence of gas-rich m ergers. Obse r vationa l evidence of 
this pro c ess is presented by Robin et all d 1 996bh. iRobin et alJ 
dl996al) . IGilmore et ail d2002l) and IWvse et alj d2006l) . After 
the merger, it is plausible that the star formation stopped for 
a while until the gas asse mbled again in the thin disk (e.g. , 
see extensiv e discussions in lFreeman & Bland-Hawthorn| [2002; 
Wyse 2002). Finally the MW can produc e thick disk f e atures 
in itself by radial mig ration processes (Roskaretal. 2008; 
Schonrich & Binnev 2009) or from disruption of massive star 
clusters (e.g.. lKroupall2002l) . 

In order to distinguish the role of these different scenar- 
ios the improved kinematics data expected from th e forthcom- 
ing as trometric Gaia satellite are fundamental (e.g., Rob in et al.l 
2012), especially in order to constrain time-e volving self- 
consistent dynamical and chemical models ( e.g., iRahimi et al 



(46.1 + 6.7)kms- 
3.4)km s" 1 . 



cr™ = (5.8 ± 5.1)kms-\ cr zz = (35.1 ± 



l20lUlLoe bman et al. 201ll lHouse et al.ll201 li M oni Bidin et al 
2012a). If the thick disk is formed from accreted stars, e.g. dur- 
ing a merger event, then no verti cal gradient, [Fe/ H] (z) is ex- 
pected for the thick disk (but see iKatz et aUl201 ll). The radia l 
mixing for a MW in isolation (e.g. Schonrich & Binnev 2009) 
is probably not an efficient mechanism to remove t he metallicity 
(or alp ha elements) radial gradients (if any, see e.g jRuchti et al.l 
1201 lal) beyond R > 9 kpc because of the ste eply decreas- 
ing p robability of radial migration suggested by iBensbv et al.l 
(1201 lh . Radial migration remains still a mechanism to be fheo- 
retically defined beyond the mere N-body n umerical experiment 
dSolwav et alj 1201 2t iMinchev et al.l 1201 2|) . while mergers can 
sensibly flatten the migration pro bability or also a nti-correlate 
it with the radial density profile (Bird et all 1201 2l) . A correla- 
tion between ||v c || and [Fe/H] would probably disf avour a migra- 
tion s cenario and slow heating mechanisms (e.g.. lNavarro et akl 
l201lh . 
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Appendix A: The whole sky symmetry 



Here we outline the procedure to pass from a radial velocity set of data to the true velocity first order moments in the case of an 
all-sky-survey coverage. This is an application of what is shown in Section lXTl that we used as reference case. In this case the matrix 
in Eqn. dTOb can be computed directly using spherical coordinates {x,y, z} = {cos b cos/, cos/? sin/, sin/?}. The generic element of 
the matrix, Mpy, 

M [2 ]l,l 



(H 



Mi 



[2] 1,6 



Mi 



[2] 1,6 



M, 



[2]6,6 ) 



can be computed as 



M 



where S 2 = [0, 2n[ x [—n/2, n/2[ and the solid angle c/Q = did (cos /?). Or, element by element: 
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(A.3) 



From its inverse we can easily obtain the elements of the velocity ellipsoid as: 

Act 2 — cr 2 — cr 2 
HU \\RR U U4> Wzz 
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(A.4) 



In general, all this is valid when the data have a spherically symmetric distribution. With real data, owing to the partial sky coverage, 
the matrix dA.3t may substantially differ from the symmetric case. To take this into account, we use the matrix ( IA.3b as a mask, 
i.e. as a constraint on the relative weight that the generic matrix element Mpj/j has with respect to another element Mpjy when 
trying to maximize the parameter q of Eqn. (1171 1. i.e. to control the coupling of the off-diagonal blocks of the operator matrix P. The 
generic matrix element to be determined requires integrals of many complex trigonometric functions, thus implying long tedious 
calculations. Fortunately, many elements are null by symmetry. For instance, because there exist only three linearly independent 
isotropic fourth-rank tensors related to the Kronecker delta tensor by J,^/ = 6u6ik, STiju = dikdji and 'Kiju = SijSu, the generic 
symmetric fourth-rank isotropic tensor can be expressed as a linear combination of these (e.g., with coefficient X s , X a , X Tr ): 
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(I + J)-- r K\ + X a (I-J) + X' r -K. 



(A.5) 



so that only the different terms can be easily singled out from Eqn. IA.5I Comparing the shape of the matrix in Eqn. (I A. 31 > with the 
true case of RAVE data from Fig. IA. II shows that in most cases nearly symmetric conditions apply. Indeed, matrices in Eqn. ([Tol l. 
M[2], Eqn. ( fTTT i. Mpj, and ( TT2b . M[4] closely resemble the case of spherical symmetry (see, e.g., Fig. lA.ll for the matrix of Eqn. 
(TTUb). For the third moment we have 
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1 2 3 4 5 6 




Fig. A.l. The inverse matrix of the Eqn. d 1 Ob for the best fit 
model. The matrix closely resembles the fully analytical case 



1 2 3 4 5 6 except for the numerical part. Different shades of colour are 

applied to visualize the symmetries. 



and from its inverse 
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and from its inverse 



T = (Trrrr, Trrrq, Trrr z , 7^00, Trr$ z , Trr zz , 7 , «000, Trqq-, r«0 ZZ , Tr zzz , T^q, T^fe, T^^, T<p zzz , T ZZZZ J - 
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Appendix B: The V-cumulants mixture distribution 

To disentangle the first and second cumulants of the thick disk alone, we need up to the fourth order cumulants of the mixture. 
The V-cumulants mixture distribution of third and fourth order are introduced here. These moments exhibit a s ymmetry along the 
azimu thal velocity vector of the mixture. Hence they differ from what is laid out in the Appendices B, C, and D of Cubarsi & Alcobel 
(12004 . Also the notation differs in order to be consistent with the notation used in our study. The third cumulants from which we 
can easily compute the two components of the normalized vector d defined in Eqn. (fTTl l. are 



KfiRR — k zzz - d^K^ + 3j 3 /c^ z - 1)d-iK^ zz - 0, 

KfiR- = <?i%j - k Rzz + d 3 {-diK RH + 2k R4>z + dydiK^ - 23^^ = 0, 
Kfe z = k RRz - d3K RRt f, + d\ (2J3KR00 — 2n Rl p z — d\d-$K,f„f„f, + d\K$^ - 0, 
K zzz = -k rrr + d\ (3k RR4 , + d\ (-3K R <p<p + Ji*^)) = 0. 



(B.l) 



Once the values of d\ and dj, are computed as SVD solution of the previous overdetermined system, the remaining cumulants can 
be calculated via 



Krr0 — a:0 ZZ — 2(i3/f00 Z + d\ 

Kr<^ = K^ z — d^K^, 
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(B.2) 
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The fourth order equation can be computed as 

K.RRRR = Krrrr + J3 (-^Kpzzz + ^3 - 4^3^00. + d^K^Aj , 

K«««z = —krzu. + 

^3 (3f«0 Z; — 3d\K,/„f, zz + f — 3kr^, z + d-iKR^fj, + 3<ii/c^ z — did^K^^ , 
Krrzz = krr zz - 2diKm<pz + d^KRR^f,+ 
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Once we have the whole set of the cumulants we define the tensor 

1 



C = 



(fthick - ""thin) - 



(B.4) 



where d = ^d. The constraining equations derived in ICubars i & Alcobe (2004) can be reduced to the following set of fourteen 
scalar relations 
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(B.5) 



that we solve in a least-squares sense with respect to the elements and di- The final step to calculate the desired results for q is 
to take the derived values for and d.2 and to work out the last constraining equation from the relations 



K,zz = K zz 



3c00O?2 + 2^2' 
Kzzzz = Kzzzz = 3c^-2(<7 2 + l)J: 



(B.6) 



Appendix C: Thick disk parameters without knowledge of photometric distances 

By testing the meth od developed in this paper on a completely synthetic catalogue created with the Padua Galaxy Model (e.g., 
IVallenari et alj|2004l and references therein) we can test our ability to recover the correct results, to refine the method and to 
improve its performance. 
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Fig. C.l. Distance distribution for the selected sample of stars. 
The number is normalized to the highest value to evidence the 
peak position at 1. The x axis shows distances divided by the 
adopted solar position (to better illustrate the parameter e for 
which the adopted approximations hold). 



Table C.l. Values of the thick disk velocity dispersion tensor without the knowledge of the distances. 





[kms" 1 ] [kms" 1 ] [kms" 1 ] 


O'RR (Tgij, <X R - 
0~$$ 0~ 0, 


56.1 ± 3.8 29.4 ±17.2 10.1 ±3.3 
46.1 ±6.7 5.8 ±5.1 
35.1 ±3.4 



Moreover we obtain here a remarkable example of convergence of the method on th e true RAVE data, w here the method is forced 
to work without the knowledge of the previously determined photometric distances by Zwitter et aT| (120101) . A fundamental selection 
criterion in order to achieve this particular result is the cut in the surface gravity of the stars. In order to avoid t he contamination by 
giant stars which can enter our sample because of their intrinsic luminosity despite their distances (see, e.g., KJ ement et al.l l2008, 
1201 ll) we plot in Figure ICTI the distance distribution of the dwarf stars selected with a cut in the surface gravity at log 10 g = 3.5. 
Moreover in the introductory consideration (Section|2| we explained how the expansion over a parameter s of the radial component 
of the Galactic rotation leads to only a weak influence of the photometric distance errors on the results. This parameter has to be 
small, of the order of e - |. As evident by plotting the distances for an averaged sample of stars of the mixture (« 38, 000 stars) the 
distribution shows an e variation well within the | = 0.12 confirming that our selection cuts are able to retain stars with distances 
within a range of r^i < 1 .Okpc without an a priori knowledge of the distance. 

After that, technically speaking, this result is achieved by simply allowing the variation of the randomly sampled distances of 
each star not only within their photometrically determined errors, but along all of the lines of sight, working only with directions 
(/, b) instead of the full parameter space of directions and distances (I, b, d). The results are presented in Table IC.fl The results are 
remarkably similar to the ones presented in Table[T]as expected from the selection criteria adopted in Eqn. (|5). 
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